### Perform spatial cross-correlation test and visualization
### Single point observation algoritm (both observation in 1 geo point)

# ref to: https://storage.googleapis.com/plos-corpus-prod/10.1371/journal.pone.0126158/1/pone.0126158.pdf?X-Goog-Algorithm=GOOG4-RSA-SHA256&X-Goog-Credential=wombat-sa%40plos-prod.iam.gserviceaccount.com%2F20210921%2Fauto%2Fstorage%2Fgoog4_request&X-Goog-Date=20210921T133943Z&X-Goog-Expires=86400&X-Goog-SignedHeaders=host&X-Goog-Signature=9de9a26c3e63c996dc8dc97816a748fbedb7cf307a2ad2ceb903a47b2d0a7b83f222e65efd1dc317af0f268db372787aa9ca9f8ddaa49b43c7e0ad1632f665bed702820e0a040b1d68bf7b8daccce3c93fa23e52f6907706c77d887e9c2535e84209200df1e17f6828a81e0b08995be273b2ce3dc36d93554606fe5910cb48f2e39704d7d95f8874f0818c233753aeb9809bd0921ce256d95b728d74f78eac3d1af67b2ff0d4c1c460e2226e6aae5272e2a433fd3d9fdbb25cf83429d625f033bf9895ffd26b9c756cd67df19e6cad4a7ca7b028c62f39be1df54b2e69ad2c94fa9afa04073369e813c04c710e5efc8cca42afff96d8fe50fc94ef07f6e6041c
# ref to: https://cran.r-project.org/web/packages/spatialEco/spatialEco.pdf
# ref to: https://gis.stackexchange.com/questions/278982/correlation-relationship-between-map-layers-in-r

# load lib
library("spatialEco")
library("sp")
library("spdep")

# load data
data <- read.csv("input/data-same.csv")
coordinates(data) <- ~x+y

# make 3 objects for different analysis purposess
data.direct <- data
data.eucl <- data

#### Example 1. Direct crossCorrelation experiment
# perform crossCorrelation test with 99 iterations 
fit.direct <- crossCorrelation(data.direct$fact1, data.direct$fact2, coords = coordinates(data.direct), k = 99)
# extract SCI statistics
data.direct$lisa <- fit.direct$SCI[, "lsci.xy"]
# plot spatial correlation
spplot(data.direct, "lisa")

### Example 2. Smoothed distance weighted matrix by Euclidian distance between points
weigths.euclidian <- spDists(data.eucl)
fit.euclid <- crossCorrelation(data.eucl$fact1, data.eucl$fact2, 
                               w = weigths.euclidian, 
                               dist.function = "none", 
                               k = 99, 
                               scale.partial = T,
                               #scale.matrix = T,
                               type = "GSCI")

data.eucl$lisa <- fit.euclid$SCI[, "gsci.xy"]
spplot(data.eucl, "lisa", main = "Spatial CC (direct point)")

# Perform cross-correlation linear weighted test
morans.plot(x = data.eucl$fact1, y = data.eucl$fact2, coords = coordinates(data.eucl), scale.morans = T)
